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One-dimensional (ID) Nonlinear Schrodinger Equaation (NLS) provides a good approximation 
to attractive Bose-Einshtein condensate (BEC) in a quasi ID cigar-shaped optical trap in certain 
regimes. ID NLS is an integrable equation that can be solved through the inverse scattering method. 
Our observation is that in many cases the parameters of the BEC correspond to the semiclassical 
(zero dispersion) limit of the focusing NLS. Hence, recent results about the strong asymptotics of the 
semiclassical limit solutions can be used to describe some interesting phenomena of the attractive 
ID BEC. In general, the semiclassical limit of the focusing NLS exibits very strong modulation 
instability. However, in the case of an analytical initial data, the NLS evolution does displays 
some ordered structure, that can describe, for example, the bright soliton phenomenon. We discuss 
some general features of the semiclassical NLS evolution and propose some new observables to the 
attractive ID BEC. 

PACS numbers: 



ID BEC as a semiclassical limit of the focusing 
NLS. It is generally accepted that the temporal evolution 
of the BEC is governed by the Gross-Pitaevskii (GP) 
equation for the condensate wave function ^(r,t), given 

by 
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where: m is the single atom mass; V ext is an external 
trapping potential that we consider to be a harmonic 
oscillator potential, and g is the coupling constant de- 
termined by the scattering length. Negative values of g 
correspond to the attractive BEC. If the characteristic 
energies of the radial excitations are much greater then 
the energy of the nonlinear term (sec estimate (5) below), 
the 3D GP equation (1) can be approximated by a ID 
GP in the longitudial (axial) direction ([!]) 
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where £ is the axial variable and 
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is the axial part of the harmonic trapping potential. Here 
loj_ and W{ are the trap frequcn- 
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cies in the radial ana the axial directions respectively. 
Taking as an example the bright solitons experiment with 
Lithium ( 7 Li, Sreckker et al., [10]), we have m w 10~ 26 kg, 
a s w — 3a rts —16- 10~ n m (here do denotes Bohr radius) 
and oj± — 2ir ■ 640Hz« 4 • 10 3 H; in the first approxima- 
tion, we put = 0. The axial wave function ^>(£, t) is 
assumed to have a shape of Gaussian and is normalized 

by 
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where N 3 • 10 5 is the total number of the atoms in the 
trap. 

Considering |*| 2 - J?-, the ID GPE (2) is applicable 
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([1]) under the condition 
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where sm is the order of magnitude of the size of the 
condensate in the axial direction. Assuming S|| ss 3 • 
10~ 4 m, ([10]), the left hand side of (5) becomes 0.16, 
which may give some justification for the use of (2). 

Equation (2) with zero external potential is a ID NLS, 
which can be integrated through the inverse scattering 
technique. We start our discussion by showing that the 
GP equation (2) with zero external potential that de- 
scribes the attractive BEC in a cigar-shaped trap can be 
rescaled to 
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where e is a small positive parameter, x, r are scaled 
space - time variables and the initial data |g(x,0)| has 
a shape of Gaussian with a typical lenght of order 1 is 
normalized by 
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| 9 (x,0)| 2 da; = 1 . 
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Equation (6) is the standard form of the focusing NLS in 
the semiclassical (zero dispersion) limit. 

Substitution of expressions for g, l± into (2) yields 
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ip(Z,t)=0. (8) 



Equation (6) can be obtained from (8) through the 
change of variables 

i/>(Z,t) = 0q(x,T), Z = Ax, t = kr, (9) 
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where the coefficients [3 1 A and k are to be determined. 
The comparison of the norming conditions (4) and (7) 
yields (i 2 — ^. Substituting (9) into (8) one gets 
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(10) 

Comparison of equations (10) and (6) yields 
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Comparison of the typical size of Gaussian distributions 
for |</>(£, 0)| and \q(x, 0)| gives 

A ~ 10~ 4 . (13) 
Using numerical value H 10~ 34 m s ks , we calculate 

ewl.6-10" 2 and k w 5 • 10~ 3 . (14) 

Equation (14) shows that the time evolution of attractive 
BEC, governed by equation (8), can be described by the 
semiclassical limit of the focusing NLS (6) with normal- 
ization (7). (Normalization of J \q\ 2 dx = n instead of 
(7), where n — 0(1), leads to the replacement of N by 
N/n in equations (10)-(12). Thus, in this case e,k has 
to be multiplied by \/n.) 

Semiclassical limit solutions to the focusing 
NLS (6). The focusing Nonlinear Schrodinger equation 
(6), where x € K and r > are space-time variables, 
is a basic model for self- focusing and self- modulation; 
it describes the evolution of the envelope of modulated 
wave in general nonlinear systems. It is also one of 
the most celebrated nonlinear intcgrable equations that 
was first integrated by Zakharov and Shabat [16], who 
produced a Lax pair for it and used the inverse scat- 
tering procedure to describe general decaying solutions 
(limui^Q q{x, 0) = 0) in terms of radiation and solitons. 

In the semiclassical limit (e — > 0) the focusing NLS 
(6) exibits modulationally unstable behavior (see Fig. 1), 
as was first shown in [6]. This is in drastic contrast to 
the case of the defocusing NLS equation ([2], [7]) in which 
the semiclassical theory shows regions of modulated peri- 
odic or quasiperiodic oscillation. These two very different 
types of behavior can be explained through modulation 
equations, which are elliptic in the focusing and hyper- 
bolic in the defocusing cases. The corresponding initial 
value problems are, therefore, ill-posed and well-posed 
respectively. As a result, a plane wave with amplitude 
modulated by A(x) and phase modulated by S(x), taken 
as an initial data 




FIG. 1: Absolute value \q(x,r,e)\ of a solution q(x,T,e) to 
the focusing NLS (6) versus x, t coordinates from [2]. Here 
A(x) = e-* 2 , S'(x) = -tanhz and e = 0.02. 



for the focusing NLS (6), is expected to break immedi- 
ately into disordered oscillations of both the amplitude 
and the phase. However, in the case of an analytic initial 
data, the NLS evolution displays some orderly structure 
instead of the disorder suggested by the modulational in- 
stability, see [2], [9] and [3]. Throughout this work, we 
will use the abbreviation NLS to mean "focusing Nonlin- 
ear Schrodinger equation" . 

Fig. 1 from [2] depicts the time evolution of a typical 
Gaussian-shaped symmetrical analytic initial data (15). 
It clearly identifies regions where different types of be- 
havior of the solution q(x, r, e) appear. These regions 
are separated by some independent of e curves in the 
x, t plane that are called breaking curves or nonlinear 
caustics. Within each region, the strong asymptotics of 
q(x, t, e) can be expressed in terms of Ricmann Theta- 
functions (see, for example, [11]). In this context, regions 
of different asymptotic behavior of q{x 1 r 1 e) corresponds 
to the different genera of the Schwarz symmetrical hyper- 
elliptic Riemann surface TZ(x,t), whose Theta-functions 
enter in the asymptotic description. In the genus zero 
case, 1Z(x,t) has only two branchpoints a — a(x,r) and 
a, so that the Riemann Theta-functions expression for 
q(x, t, e) is replaced by 



q(x,r,e) = A(x,T)e lS(x ' T)/£ , 



where 



a(x, t) 



-S x (x,t) +iA(x,r). 



(16) 



(17) 



q(x,0,e) = A{x)e lS( - x)/e 



(15) 



The genus zero region is the first region adjacent to the 
axis r — 0, where (genus zero) solution (16) has the form 
of a high frequency (0(l/e)) modulated wave with the 
slowly varying amplitude A(x,t) and phase S(x,t). (We 
remark that A(x,0) = A(x), S(x,0) = S(x).) 

In some very few special cases, for example, when 
A(x) = sech x and S'(x) — — ^tanhx, fi > 0, the scat- 
tering data for (6) can be calculated explicitly. Then 
the modulated amplitude and phase of (16), (17) can be 
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obtained from the system of trancendental equations for 
a{x, t) — a(x, t) + ib(x, r): 



y/(a - T) 2 + b 2 + ^(a + T) 2 + b 2 = n + 4r6 2 
a-T+ ^{a-T) 2 + b 2 } \a + T+ ^{a + T) 2 + b 2 

(18) 

where T = \J ^ — 1. In the particular case /i = 2 (the 
borderline value of /i between the pure radiational case 
fi > 2 and radiation with solitons case fi < 2), introduc- 
ing the implicit time u — u(x, r) at each point x € K 

i (u— :r)[sinh 2u — (u— x)] i < • i- • < 

by t = - 8sinh 2 u > one can 0D t am an explicit 

solution ([11]) 



2 sinh 2 u 
sinh 2it — (u — x) ' 
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sinh 2u — (u — x) 



for A(a;,r) = 6 and S"(x,t) = —2a. Similar expressions 
are available for the case fi = 0. 

Notice that the amplitude A(x, t) of the solution on 
Fig. 1 at first contracts (accumulates) towards the point 
of maximum (x = 0) of \q(x, 0)| and then suddenly bursts 
into rapid (order 1/e) and violent oscillations in ampli- 
tude (transition to genus two regime). This is the typi- 
cal behavior ([12]) for an analytic one-hump initial data 
(provided that S'(x) does not decrease too fast) The very 
first point of this transition, which is the tip-point of the 
first breaking curve (see Fig. 1) , is called a point of gradi- 
ent catastrophe, or elliptic umbilical singularity ([4]). At 
the point (xq,tq) of gradient catastrophe the semiclassi- 
cal solution (16) of (6) losses its smoothness ([15]), i.e., 
a x (x ,T ) = oo (either A x (x,tq) or S xx (x,t ) or both 
become infinite). 

The Theta-function expression for higher genus solu- 
tions is somewhat cumbersome for this paper (see, for 
example, [11]). This expression gives an O(e) approxi- 
mation for the solution of (6), (15) in the corresponding 
region. For a fixed time snapshot in the genus two re- 
gion (Fig. 1, around r = 1), the graph of \q(x,r,e)\ can 
be identified with bright solitons, that were experimen- 
tally observed, for example, in [10]. If the BEC with a 
Gaussian shaped initial data is governed by ID NLS (8), 
the region filled with solitons is spreading off in the axial 
direction. According to Fig. 1, the onset of soliton-likc 
(genus two) behavior in the semiclassical limit happens 
at r = \ or t = kr « 2.5 • 10~ 3 s, which does not contra- 
dict observations of [10], were bright solitons were first 
observed at t = 5ms. Since the "effective" axial size of 
the condensate sn shrinks considerably near the time of 
gradient catastrophe (see Fig. 1), condition (5) may be 
violated during this period. This is consistent with the 
fact that the total number of atoms observed in the soli- 
ton regime in the experiment of [10] is less than 20% of 
the number of atoms N m 3 • 10 5 at the beginning of the 
experiment. (The NLS evolution preserves the L 2 norm 
of the solution, i.e., evolution governed by equation (8) 
would preserve the total number of atoms N.) Moreover, 



Momentum distribution as t increases in the 
Genus region 




♦ ♦ 



+ t=.125 
■ t=.100 

Xt=.D5Q 

«t=D 



relative momentjrr iPj 



FIG. 2: Evolution of the Fourier transform of the initial data 
(15) with A(x) — sech x and ^'(a;) = — 2tanha; in the limit 
e — > from r = to the time of gradient catastrophe r = 
0.125. 



the evolution of the Fourier transform of q{x, t, e), which 
can be calculated explicitly form (16) through the sta- 
tionary phase method, shows that the portion of atoms in 
the condensate with high axial momentum significantly 
increases (see Fig. 2) as the point of gradient catastrophe 
is approached. 

Note that breaking curves (boundaries between the re- 
gions of different genera) depend on both the amplitude 
and the phase of the initial condition q(x,0,e), but do 
not depend on e. Within a region of genus 2n, n > 0, 
the semiclassical limit solution can be viewed as mod- 
ulated 2n- phase nonlinear wave, with 2n + 1 complex 
(4n + 2 real) wave parameters, which slowly vary (in 
x, r) in the region. (Fig. 1 depicts consequitive re- 
gions with n — 0,2,4). Complex wave parameters can 
be interpreted as a set of branchpoints of the Schwarz- 
symmetrical hyperelliptic surface lZ(x, r), whose evolu- 
tion in the x, r plane is defined by modulation (Whitham) 
equations. (Here we want to mention that regions of 
genus four (see Fig. 1 after r = 1.5) and of higher gen- 
era are associated with the initial data (15) that support 
solitons, see [11]; in the semiclassical limit solutions the 
number of solitons has order 0{l/e).) 

Calculation of semiclassical solutions to (6). 
Equation (6), as an integrable NLS, can be solved by 
inverse scattering technique. However, the semiclassical 
limit solutions require the semiclassical limit of the scat- 
tering transform. Let z be a point on the curve E in the 
upper halfplane, that is defined parametrically by the an- 
alytic initial data (15) as a(x) = — \S'{x) + iA{x), x e M. 
(Here A(x) and S'(x) have sufficient decay to zero or to 
some finite values ±/z at ±oo respectively.) Assuming for 
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simplicity that a(x) is invertible, the semiclassical scat- 
tering data limit fo(z), z e X, is defined ([15]) through a 
generalized Abel integral transform as 



fo(z) 



z — ii + + \J \z — u)(z — u) x'(u)du 



+ {z-n+)x{z), (20) 

where x(a) is inverse to a(x) and the integral is taken 
along S. The analytic extension of fo{z) from E to 1 
(which can have logarithmic branchcuts) has a meaning 
of the leading order term of \ie In r (z, e) as e — > 0, where 
ro(z,e), z € E, is the reflection coefficient of (15). Once 
/o(z) is known, the complex wave parameters arc defined 
through the modulation equations. In particular, in the 
genus zero region, the modulation equation for a(x, r) is 
given by a system of two real equations ([11]) 
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where f(z) 

sj(z — a){z — a). It defines q(x,r,e) through (16)-(17). 
(Here fo(z) is Schwarz symmetrically extended into the 
lower halfplanc; typically, ^sfo(z) has a jump along E.) 
Define function h(z) = h(z;x,r) as 

i r /(O 



Z7T 
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(22) 
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where j m is a Shwarz-symmetrical contour connecting 
a and a, and such that 7 TO U E = Because of the 
analyticity of f{z), a particular shape of j m is not impor- 
tant. However, it is possible to fix j m by the condition 
3/i(z) = on 7 m . According to the Dcift-Zhou nonlin- 
ear steepest descent method, the genus zero anzatz (16) 
approximates the actual solution of the NLS (6) with the 
reflection coefficient r Q (z,e) — e~^°( z ) if ([11]) 

3/i(z; x, t) < on both sides of 7+; 

Qh(z;x,r) > on 7+, (23) 

where 7+ is a contour in the upper halfplane C + connect- 



ing a and and 7„ 



7r, 



U C + . We have a freedom 



to deform the contour 7+ so that the inequalities (23) 
are satisfied along it. The first breaking curve consists 
of points (x, t) where at least one of the inequalities (23) 
turns into equality at some z . Thus, equation for the 
first breaking curve can be written as a system of three 
real equations for z eC and (x, r) € E 2 

?Sh(z ; x, t) = and h z (z ; x,t) = . (24) 

For the initial data (15) with when A(x) = sech x and 
5'(x) = — tanhx, the expression 

Va 2 + b 2 R(z) - a(z - a) + b 2 



h(z) = zln - 
+ (1-*) 



- 2t(z - a)R(z) 
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ln[R(z) -{z- a)} 
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was found in [11]. Modulation equations, as well as ex- 
pressions for h(z;x,r) for higher genus regions, can be 
written in the explicit determinantal form (see [13], [14]), 
however, since these expressions are somewhat involved, 
they will not be given in this paper. 

Suggestions and conclusions. Semiclassical limit of 
the focusing ID NLS (6) provides a new, mathematically 
rigorous tool to study modulationally unstable evolution 
of the attractive ID BEC. Evolution of the BEC with a 
one-hump initial data that is governed by the NLS (8) 
is expected to show two or more qualitatively different 
regimes (regions of different genera) within 0{k) (see eq. 
(12)) time interval. Bright soliton experiment of [10] is 
an example of a typical higher genera region behavior 
(apparently, for the genus 2 region), attained within the 
time period of 5ms, where k 4ms. In general, any 
macroscopic characteristic of the evolving condensate can 
be suggested as an observable. That include (see Fig. 1): 

• space-time location of the breaking curves (exact 
location of a breaking curve is given by (24) ) ; 

• slowly modulated amplitude A(x,t) in the genus 
zero (the exact value of A(x,t) is given through 
the modulation equation (21)); 

• the upper and lower envelopes of the high-frequency 
amplitude oscillations in the genus two region (the 
envelopes are defined through the Ricmann-thcta 
functions) . 

Calculation of the observables, mentioned above, is 
based on the semiclassical limit of the scattering data 
fo{z), which, in its turn, can be obtained from q(x, 0,e) 
through (20), i.e., through the initial amplitude A{x) and 
the phase S(x). However, accurate measurment of the 
initial phase is often a difficult task. We can turn the 
question around and ask whether the phase S(x) can 
be somehow reconstructed from A(x) and some observ- 
ables. Continuation of this line of argument leads to the 
question of designing some NLS-data, initial or scatter- 
ing, whose evolution will have certain desired properties 
and/or fit within some required parameters. Generally 
speaking, formulae (20)-(24) are valid for a large class 
of analytic initial data, including, for example, cases of 
the multi-hump initial data A{x), and experimental data 
about evolution of the BEC with multi-hump initial den- 
sity might be interesting. 

ID NLS approximation of the evolving attarctive BEC 
in a cigar-shaped trap may be valid during certain in- 
tervals of the total period of observation and not valid 
during the others. The atoms lost in the experiment of 
[10] (probably, near the point of the gradient catastro- 
phe), seem to indicate that this is an example of such 
situation. Perhaps, some other model is needed to trace 
the evolution of the BEC through the point of the gra- 
dient catastrophe into the soliton regime, where the ID 
NLS approximation will be working again. 
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